Einführung. Erklärung. Schlussfolgerungen mit Links zu den relevanten Teilen.
1.) Stimmensplitting
2.) P-Werte
3.) Signifikanz
4.) Wahlbeteiligung Teil 1 und Wahlbeteiligung Teil 2

Unter jedem Titel steht "Hoch." Wenn man da drauf klickt, gelangt man zurück zum Inhaltsverzeichnis.

Inhaltsverzeichnis¶

1. Daten aufbereiten
   • Wahldaten 2021
   • Einkommensteuer-Verteilung 2020
   • Ausländer und Gesamtbevölkerung
   • Arbeitslose und Hartz IV
   • Alles Zusammen
2. Daten analysieren
   • Einfache Korrelationen und Plots
     a) P-Werte
     b) Multiple Regression)
     c) Ausreißer?)
     d) Variance Inflation Factor)
     e) ANOVA)
   • Faktorenanalyse
     a) Sind die Daten für Faktorenanalyse geeignet?)
     b) Wie viel Faktoren?)
     c) Regression mit Factor Scores)
     d) Signifikanz
     e) Auswertung)
   • LASSO und Validierung
     a) Das LASSO-Modell anpassen)
     b) Interpretation des LASSO-Modells)
     c) Validieren der Modelle)
     d) Verstehen vs. Voraussagen)
   • Mit Kanonen auf Spatzen: Random Forests und Boosting
3. Was könnte man sonst noch machen?

import pandas as pd
import numpy as np

from factor_analyzer import FactorAnalyzer, calculate_kmo, calculate_bartlett_sphericity

import matplotlib.pyplot as plt
from matplotlib.pyplot import subplots

from sklearn.ensemble import RandomForestRegressor
from sklearn.linear_model import LassoCV, LinearRegression
from sklearn.metrics import mean_squared_error
from sklearn.model_selection import train_test_split
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import StandardScaler

from scipy.stats import pearsonr

import statsmodels.api as sm
import statsmodels.formula.api as smf
from statsmodels.stats.anova import anova_lm
from statsmodels.stats.outliers_influence import variance_inflation_factor as VIF

import xgboost as xgb

Daten aufbereiten ¶

Hoch

In diesem Teil geht es nur darum, wie die Daten verarbeitet werden mussten, um nutzbar zu sein. Er kann ohne Verständnisverlust übersprungen werden.

Alle Daten sind vom Landesamt für Statistik, aber sie sind nicht unbedingt im nutzerfreundlichsten Format, weshalb ich alle Daten ein wenig bearbeiten musste und sie nicht einfach direkt ins Notebook laden kann. Mehr zu den einzelnen Verarbeitungen weiter unten, aber zusammengefasst:

1.) Die Wahldaten-Datei enthält immer alle Daten für alle Regionen, auch wenn man versucht, nur die Stadtteildaten herunterzuladen. Daher müssen einige Zeilen gelöscht werden. Die erste Zeile enthält nie die Spalten-Titel, sondern die Überschrift der Tabelle. Daher müssen entweder ein paar Zeilen übersprungen oder alle Spalten umbenannt werden.

2.) Aus der Einkommensteuer-Verteilung für das Land Bremen kann man ablesen, wo die Grenzen für das unterste und oberste Fünftel liegen. Für jeden einzelnen Stadtteil rechne ich daher aus, wie viele jeweils unterhalb und oberhalb dieser Grenze leben. Das statistische Landesamt stellt diese Daten leider nicht direkt bereit.

3.) Auch den Anteil an der Gesamtbevölkerung in den Stadtteilen, den Ausländer darstellen, bietet das Landesamt, glaube ich, nicht direkt an. Ich errechne ihn daher selbst.

4.) Die Arbeitslosenquote wird wohl nur in Quartalen angeboten. Ich errechne daher die durchschnittliche Quote im Jahr vor der Bundestagswahl 2021. Dasselbe gilt für die absoluten SGB-II Leistungen, die nur in Quartalen verfügbar sind. Außerdem errechne ich hier auch die Leistungen pro Kopf, da die absolute Zahl v. a. von der Bevölkerungszahl abhängt.

Wahldaten 2021 ¶

Hoch

Erst- und Zweitstimmen sind auf unterschiedlichen Blättern in der Excel-Datei, können aber fast gleich behandelt werden. Die Zahlen zur Briefwahl sind nur auf dem ersten Blatt und müssen nur davon gelöscht werden. Außerdem ignoriere ich kleine Parteien. Nur SPD, CDU, Grüne, Linke, FDP und AfD bleiben im DataFrame^†. Die Regionen Seehausen, Strom, Borgfeld und Blockland sind zwischen den Stadtteilen und müssen gesondert entfernt werden.

† Ein "DataFrame" ist eine mächtigere Version einer Excel-Tabelle. Verschiedene Arten von Daten (Zahlen, Text, Kalenderdaten, ...) können in seinen Zeilen und Spalten gespeichert werden und man kann diese Daten leichter analysieren, filtern oder sonst verarbeiten

erststimmen_2021 = pd.read_excel('Wahldaten 2021.xlsx', sheet_name=0) 

# Fixing the wrongly imported column names.
# Because of the layout, if you don't set that NaN value at [0, 1] to 'region' (or any other not NaN), both the first and the second column get dropped in the next step.
erststimmen_2021.iloc[0,1] = 'region'
erststimmen_2021.columns = erststimmen_2021.iloc[0]

# Dropping Codes, Briefwahl (4 and 5), stärkste Partei (6 and 7), and small parties.
drop_columns = [0, 4, 5, 6, 7] + list(range(20, erststimmen_2021.shape[1]))
erststimmen_2021 = erststimmen_2021.drop(erststimmen_2021.columns[drop_columns], axis=1)

# Isolating the Stadtteile
start_index = erststimmen_2021[erststimmen_2021['region'] == 'Mitte (Bremen)'].index[0]
df_filtered_partial = erststimmen_2021.loc[start_index + 1:]
end_index = df_filtered_partial[df_filtered_partial['region'] == 'Blumenthal'].index[0]

df_filtered = erststimmen_2021.loc[start_index:end_index]

# Dropping rows of non-Stadtteile in between
df_filtered = df_filtered[~df_filtered['region'].isin(['Seehausen', 'Strom', 'Borgfeld', 'Blockland'])]

# Adding in Land Bremen, Stadt Bremen, Stadt Bremerhaven
df_add = erststimmen_2021.iloc[[161, 159, 160]]
erststimmen_2021_sortiert = pd.concat([df_add, df_filtered], ignore_index=True)

# Making the names consistent with my other tables
cleaned_erststimmen_2021 = erststimmen_2021_sortiert
cleaned_erststimmen_2021['region'] = erststimmen_2021_sortiert['region'].replace({'Mitte (Bremen)':'Mitte', 'Oestliche Vorstadt':'Östliche Vorstadt'})

zweitstimmen_2021 = pd.read_excel('Wahldaten 2021.xlsx', sheet_name=1)
zweitstimmen_2021.iloc[0,1] = 'region'
zweitstimmen_2021.columns = zweitstimmen_2021.iloc[0]

# ... Different numbers but same thing
drop_columns = [0, 2, 3] + list(range(16, zweitstimmen_2021.shape[1]))
zweitstimmen_2021 = zweitstimmen_2021.drop(zweitstimmen_2021.columns[drop_columns], axis=1)

# ...
start_index = zweitstimmen_2021[zweitstimmen_2021['region'] == 'Mitte (Bremen)'].index[0]
df_filtered_partial = zweitstimmen_2021.loc[start_index + 1:]
end_index = df_filtered_partial[df_filtered_partial['region'] == 'Blumenthal'].index[0]

df_filtered = zweitstimmen_2021.loc[start_index:end_index]

# ...
df_filtered = df_filtered[~df_filtered['region'].isin(['Seehausen', 'Strom', 'Borgfeld', 'Blockland'])]

# ...
df_add = zweitstimmen_2021.iloc[[161, 159, 160]]
zweitstimmen_2021_sortiert = pd.concat([df_add, df_filtered], ignore_index=True)

# ...
cleaned_zweitstimmen_2021 = zweitstimmen_2021_sortiert
cleaned_zweitstimmen_2021['region'] = zweitstimmen_2021_sortiert['region'].replace({'Mitte (Bremen)':'Mitte', 'Oestliche Vorstadt':'Östliche Vorstadt'})

cleaned_wahl_2021 = pd.merge(cleaned_erststimmen_2021, cleaned_zweitstimmen_2021, on='region')

Einkommensteuer-Verteilung 2020 ¶

Hoch

Das statistische Landesamt benutzt das veraltete ISO-8859-1 encoding für seine Symbole statt des neueren UTF-8, das mehr Symbole repräsentieren kann und grundsätzlich von Pandas – der Library, die die DataFrame-Funktionalität bietet – erwartet wird. Außerdem kann Pandas nicht mit den Umlauten umgehen, was aber kein Fehler des Landesamtes ist.

df_einkommen = pd.read_csv('Einkommensteuer-Verteilung 2020.csv', 
                           encoding='ISO-8859-1', 
                           skiprows=2, 
                           on_bad_lines='skip', 
                           sep=';')

# Last rows are just copyright. First column are incomprehensible codes for each region. Third column is the year; 2020 for all of them.
df_einkommen= df_einkommen.drop(df_einkommen.index[242:]).drop(df_einkommen.columns[0], axis=1).drop(df_einkommen.columns[2], axis=1).reset_index(drop=True)

# The column names were a little unwieldy or contained umlauts.
df_einkommen.rename(columns={
    df_einkommen.columns[0]: 'region',
    df_einkommen.columns[1]: 'range_euro',
    df_einkommen.columns[2]: 'num_steuerpflichtige'}, inplace=True)

# More umlauts in the regions
df_einkommen.iloc[:, 0] = df_einkommen.iloc[:, 0].replace({'Häfen': 'Häfen',
                                                           'Ã\x96stliche': 'Östliche',
                                                           'Gröpelingen': 'Gröpelingen'}, regex=True)

df_einkommen.head()

Wie diese Daten direkt zeigen, haben im Land Bremen jeweils etwa 20% der Bevölkerung unter 10.000€ und über 50.000€ zu versteuerndes Einkommen. Als nächstes produziere ich daher eine Tabelle, die für jeden Stadtteil den Anteil der Bevölkerung im untersten und obersten Einkommens-Quintil anzeigt (tatsächlich sind es natürlich Einkommensteuer-Quintile).

brackets_lo = ['0-5000', '5000-10000']
brackets_hi = ['50000-125000', '125000 oder mehr']

results_arm = []
results_reich = []

counties = df_einkommen['region'].unique()

for county in counties:
    
    county_data = df_einkommen[df_einkommen['region'] == county].copy()
    
    # num_steuerpflichtige is a bunch of strings.
    # The errors='coerce' forces the missing values in "Häfen" to be NaN
    county_data['num_steuerpflichtige'] = pd.to_numeric(county_data['num_steuerpflichtige'], errors='coerce')
    
    lo_bracket_taxpayers = county_data[county_data['range_euro'].isin(brackets_lo)]['num_steuerpflichtige'].sum()
    hi_bracket_taxpayers = county_data[county_data['range_euro'].isin(brackets_hi)]['num_steuerpflichtige'].sum()
    
    # Calculate the proportion
    total_taxpayers = county_data.iloc[-1]['num_steuerpflichtige']
    
    proportion_lo = lo_bracket_taxpayers / total_taxpayers
    proportion_hi = hi_bracket_taxpayers / total_taxpayers
    
    results_arm.append({'region': county, 'anteil_erstes_quint': proportion_lo})
    results_reich.append({'region': county, 'anteil_letztes_quint': proportion_hi})

df_quint_arm = pd.DataFrame(results_arm)
df_quint_reich = pd.DataFrame(results_reich)

df_anteil_einkommen_unsortiert = pd.merge(df_quint_arm, df_quint_reich, on='region')

Dann ist noch die Zeile "Stadt Bremerhaven" ganz am Ende. Die nächste Zelle sortiert das um und entfernt das "(Stadtteil)" Addendum aus allen Zeilen.

last_row = df_anteil_einkommen_unsortiert.iloc[-1]

# Drop the last row
df = df_anteil_einkommen_unsortiert.iloc[:-1]

# Insert the last row as the second row
df = pd.concat([df_anteil_einkommen_unsortiert.iloc[:2], pd.DataFrame([last_row], columns=df.columns), df.iloc[2:]], ignore_index=True)

# Cleaning up the names. If expand=True, extract() returns a DataFrame which has not .str.strip() method.
df['region'] = df['region'].str.extract(r'([^(]+)', expand=False).str.strip()


df_anteil_einkommen = df

Da andere Daten nur für Bremen Stadt verfügbar sind, produziere ich auch eine Tabelle ohne Bremerhaven und Land Bremen.

cleaned_einkommen = df_anteil_einkommen.drop([0, 2]).reset_index(drop=True)

Ausländer und Gesamtbevölkerung ¶

Hoch

Hier musste nur der Anteil der ausländischen Bevölkerung aus den absoluten Zahlen errechnet werden.

df_ausländer = pd.read_csv('Anzahl Ausländer (Ende 2020).csv', 
                           encoding='ISO-8859-1', 
                           skiprows=1, 
                           on_bad_lines='skip', 
                           sep=';')
df_ausländer = df_ausländer.drop(['Gebietsschlüssel', 'Zeit', 'Migrationsstatus 1)', 'Bevölkerung 2).1',	'Bevölkerung 2).2'], axis=1)[3:23]
df_ausländer.columns = ['region', 'num_ausländer']
df_ausländer.reset_index(drop=True);

df_bevölkerung = pd.read_csv('Bevölkerungsanzahl 2020.csv', 
                           encoding='ISO-8859-1', 
                           skiprows=1, 
                           on_bad_lines='skip', 
                           sep=';')
df_bevölkerung = df_bevölkerung.drop(df_bevölkerung.columns[[0, 2, 3]], axis=1)[1:21]
df_bevölkerung.columns = ['region', 'num_bevölkerung']
df_bevölkerung.reset_index(drop=True)

merge = pd.merge(df_bevölkerung, df_ausländer, on='region')
merge['region'] = merge['region'].str.extract(r'([^(]+)', expand=False).str.strip()

# Convert the columns to numeric types
merge['num_bevölkerung'] = pd.to_numeric(merge['num_bevölkerung'], errors='coerce')
merge['num_ausländer'] = pd.to_numeric(merge['num_ausländer'], errors='coerce')


merge['anteil_ausländer'] = merge['num_ausländer'] / merge['num_bevölkerung']

cleaned_bevölkerung = merge
cleaned_bevölkerung

Arbeitslose und Hartz IV ¶

Hoch

Die durchschnittliche Arbeitlosenquote pro Stadtteil im Jahr vor der Bundestagswahl 2021. Das wird nicht direkt vom statistischen Landesamt zur Verfügung gestellt, sondern muss per Hand errechnet werden. Da Pandas von einem Amerikaner grundsätzlich für Amerikaner produziert wurde, kann es nicht mit der deutschen Konvention umgehen, Dezimalstellen mit Kommas statt mit Punkten abzutrennen (1,0 statt 1.0).

df_arbeit = pd.read_csv('Arbeitslosenquote (Sept. 2020 bis Juni 2023).csv', 
                           encoding='ISO-8859-1', 
                           skiprows=1, 
                           on_bad_lines='skip', 
                           sep=';')

stadtteile = df_arbeit['Gebietseinheit'].unique()
quartale = ['2021-06-30', '2021-03-31', '2020-12-31', '2020-09-30']
results_quoten = []

for stadtteil in stadtteile:
    data = df_arbeit[df_arbeit['Gebietseinheit'] == stadtteil].copy()

    data['Arbeitslosenziffer (in %)'] = data['Arbeitslosenziffer (in %)'].str.replace(',', '.')
    data['Arbeitslosenziffer (in %)'] = pd.to_numeric(data['Arbeitslosenziffer (in %)'], errors='coerce')

    avg_quote = data[data['Zeit'].isin(quartale)]['Arbeitslosenziffer (in %)'].sum() / 4.0

    results_quoten.append({'region':stadtteil, 'arbeitslosenquote_vor_wahl':avg_quote})

cleaned_arbeitslosenquote = pd.DataFrame(results_quoten)[1:21]
cleaned_arbeitslosenquote.reset_index(drop=True);

df_sgbII = pd.read_csv('SGB-II Leistungen, Sept 2020 bis Juni 2023.csv', 
                           encoding='ISO-8859-1', 
                           skiprows=1, 
                           on_bad_lines='skip', 
                           sep=';')

stadtteile = df_sgbII['Gebietseinheit'].unique()
quartale = ['2021-06-30', '2021-03-31', '2020-12-31', '2020-09-30']
results_euro = []

for stadtteil in stadtteile:
    data = df_sgbII[df_sgbII['Gebietseinheit'] == stadtteil].copy()

    data['Grundsicherung für Arbeitsuchende (SGB II) 1).2'] = pd.to_numeric(data['Grundsicherung für Arbeitsuchende (SGB II) 1).2'], errors='coerce')

    avg_euro = data[data['Zeit'].isin(quartale)]['Grundsicherung für Arbeitsuchende (SGB II) 1).2'].sum() / 4.0

    results_euro.append({'region':stadtteil, 'sgbII_leistungen_vor_wahl':avg_euro})

cleaned_sgbII = pd.DataFrame(results_euro)[1:]
cleaned_sgbII.reset_index(drop=True)

merge = pd.merge(cleaned_arbeitslosenquote, cleaned_sgbII, on='region')
merge['region'] = merge['region'].str.extract(r'([^(]+)', expand=False).str.strip()
cleaned_soziales = merge
cleaned_soziales

Alles Zusammen ¶

Hoch

Um später nicht mit einer Handvoll verschiedener DataFrames arbeiten zu müssen, Kombiniere ich hier alle sozialen Daten in einen DataFrame und auch alle Daten insgesamt in den Wahldaten_2021 DataFrame.

Stadtteile = pd.merge(cleaned_bevölkerung, cleaned_einkommen, on='region')
Stadtteile = pd.merge(Stadtteile, cleaned_soziales, on='region')

# sgbII tracks population. Creating a per capita figure.
Stadtteile['sgbII_pro_kopf'] = Stadtteile['sgbII_leistungen_vor_wahl'] / Stadtteile['num_bevölkerung']


Stadtteile = Stadtteile.set_index('region').apply(pd.to_numeric, errors='coerce')
Stadtteile

Wahldaten_2021 = pd.merge(Stadtteile, cleaned_wahl_2021, on='region')
Wahldaten_2021 = Wahldaten_2021.set_index('region').apply(pd.to_numeric, errors='coerce')
Wahldaten_2021

Daten analysieren ¶

Hoch

Ich analysiere die Daten in aufsteigender Komplexität.
    Der einfachste Schritt ist, sich erstmal Plots verschiedener Variablen gegeneinander und Korrelationen untereinander anzugucken. Das mache ich für die sozialen Daten (Bevölkerung, Sozialleistungen etc. statt Wahlergebnisse für die unterschiedlichen Parteien) und es stellt sich heraus, dass sie alle stark mit dem Erststimmenanteil der SPD korrelieren; die Korrelationen in den Daten werden höher sein, als in der Realität, da es sich um „ecological correlations“ handelt, also Korrelationen zwischen Gruppendurchschnitten.
    Es stellt sich aber auch heraus, dass alle Variablen stark untereinander korrelieren (Multikollinearität). Die „Power“ einer normalen Regression wird dadurch verringert, also die Wahrscheinlichkeit, signifikante Koeffizienten aufzuspüren; das macht die Interpretation der tatsächlich erhaltenen Koeffizienten schwierig führt aber auf natürliche Weise zu den nächsten beiden Schritten.

Da alle sozialen Variablen in meinen Daten mehr oder weniger stark miteinander korrelieren, könnte man annehmen, dass sie alle von unterliegenden, latenten Faktoren beeinflusst werden. Klassische Beispiele für so einen Zusammenhang sind psychometrische Konstrukte wie IQ und die „Big Five:“ Alle kognitiv anspruchsvollen Aufgaben – Leseverständnis, Kopfrechnen, Matrizen vervollständigen, etc. – korrelieren miteinander, weshalb Psychologen davon ausgehen, dass es einen unterliegenden allgemeinen Faktor der Intelligenz gibt, der alle speziellen kognitiven Fähigkeiten beeinflusst. Zweitens werden Menschen, die als ordentlich bezeichnet werden, auch häufiger als hart arbeitend bezeichnet. Im Allgemeinen korrelieren Beschreibungen verschiedener Persönlichkeitsmerkmale untereinander, allerdings ist ein Modell mit fünf Faktoren besser als mit einem, im Unterschied zur Intelligenz.
    Eine Faktorenanalyse der sozialen Variablen ergibt, dass es wahrscheinlich einen statt mehrere unterliegende Faktor gibt. Dieser Faktor erklärt etwa die Hälfte der Unterschiede in Erst- und Zweitstimmen unter den Stadtteilen für die SPD und AfD und wenigstens ein Viertel für die FDP und die Grünen, allerdings mit anderem Vorzeichen. Für eine nuanciertere Erklärung, was dieses Ergebnis bedeutet, s. u. unter „Faktorenanalyse.“

Eine andere Herangehensweise ist, einfach die wichtigsten Variablen auszulesen („feature selection“). Da die Variablen miteinander korrelieren, gehen kaum Information verloren, wenn man nicht alle Variablen mit ins Modell nimmt. Es gibt mehrere Möglichkeiten, die „wichtigsten“ Variablen auszusuchen, ich entscheide mich für LASSO-Regression (Least Absolute Shrinkage and Selection Operator), was komplizierter klingt, als es ist: Letztendlich wird das Modell dafür bestraft, größere Koeffizienten zu haben; einige Koeffizienten werden dadurch auf genau Null gebracht und nur die unabhängigen Variablen, die den größten Einfluss auf die abhängige Variable haben, bleiben im Modell.
    Es zeigt sich, dass der Wettbewerb um Zweitstimmen wichtiger (oder wenigstens intensiver) ist, als der um Erststimmen … nur nicht bei der Linken.

Als letztes benutze ich die stärksten machine learning Methoden neben Deep Learning: Random Forests und Boosting. Hier gibt es keine tieferen Einsichten, da beide Methoden bei diesen Daten tatsächlich schlechter funktionieren als einfache lineare Modelle. Das spricht dafür, dass es in der Realität keine komplizierten nicht-linearen Effekte oder Interaktionen zwischen den Variablen gibt.

Einfache Korrelationen und Plots ¶

Hoch

Die nächsten drei Plots zeigen die Abhängigkeit des Erststimmenanteils der SPD in der Bundestagswahl 2021 im Lande Bremen von der Wahlbeteiligung, dem Anteil der ausländischen Bevölkerung und – offensichtlich nicht kausal – dem Zweitstimmenanteil.

daten = Wahldaten_2021.drop(['Häfen', 'Stadt Bremen'])

x1 = 'Wahlbeteiligung (%) 2021'
y = 'SPD - Erststimmenanteil (%)  2021'
correlation, p_value = pearsonr(daten[x1], daten[y])
print(f'Correlation: {correlation:.3f}')
print(f'p-value: {p_value:.3f}')

plt.figure(figsize=(10, 6))
plt.scatter(daten[x1], daten[y])

regions = daten.index
for i, region in enumerate(regions):
    plt.text(daten[x1][i], daten[y][i], region, fontsize=8, ha='right')

plt.ylim(0, daten[y].max() * 1.1)

plt.xlabel(x1)
plt.ylabel(y)


plt.show()

Correlation: -0.905
p-value: 0.000

x2 = 'anteil_ausländer'
y = 'SPD - Erststimmenanteil (%)  2021'
correlation, p_value = pearsonr(daten[x2], daten[y])
print(f'Correlation: {correlation:.3f}')
print(f'p-value: {p_value:.3f}')

plt.figure(figsize=(10, 6))
plt.scatter(daten[x2], daten[y])

regions = daten.index
for i, region in enumerate(regions):
    plt.text(daten[x2][i], daten[y][i], region, fontsize=8, ha='right')

plt.ylim(0, daten[y].max() * 1.1)

plt.xlabel(x2)
plt.ylabel(y)


plt.show()

Correlation: 0.695
p-value: 0.001

Der letzte Plot (unten) zeigt eine Korrelation, wie man sie aus physikalischen Experimenten statt sozialwissenschaftlichen Daten erwartet. Anscheinend machen die Bremer keinen Unterschied zwischen Erst- und Zweitstimme, wenn sie SPD-Wähler sind; es gibt kaum Stimmensplitting.

x3 = 'SPD - Zweitstimmenanteil (%)  2021'
y = 'SPD - Erststimmenanteil (%)  2021'
correlation, p_value = pearsonr(daten[x3], daten[y])
print(f'Correlation: {correlation:.3f}')
print(f'p-value: {p_value:.3f}')

plt.figure(figsize=(10, 6))
plt.scatter(daten[x3], daten[y])

regions = daten.index
for i, region in enumerate(regions):
    plt.text(daten[x3][i], daten[y][i], region, fontsize=8, ha='right')

plt.ylim(0, daten[y].max() * 1.1)

plt.xlabel(x3)
plt.ylabel(y)


plt.show()

Correlation: 0.988
p-value: 0.000

a) P-Werte ¶

Die Korrelation und der dazugehörige P-Wert stehen jeweils über den Plots. Der P-Wert gibt in diesem Kontext an, wie Wahrscheinlich es ist, eine Korrelation zu kriegen, die wenigstens so stark ist, obwohl es tatsächlich keinen Zusammenhang gibt – dass ich also per Zufall diese Zahl erhalten habe.
    Wie man sieht kann ich ihn errechnen, obwohl er in diesem Zusammenhang bedeutungslos ist: Ich schaue mir deskriptive Statistiken der Gesamtbevölkerung an; keine Inferenz von einer Stichprobe auf die Gesamtbevölkerung findet statt, die einen Signifikanztest bräuchte. Das Durchschnittseinkommen z. B. braucht auch keinen P-Wert. Man kann ihn einfach errechnen. Wenn ich vom Durchschnitt einer Stichprobe auf den Durchschnitt der Gesamtbevölkerung schließen will, dann brauche ich einen P-Wert.

Die P-werte könnten Sinn ergeben, wenn ich die 2021 Wahldaten als Stichprobe aus einer Superpopulation über die Zeit hinweg betrachte. Die Wahldaten 2017, 2021 und 2025 wären z. B. Stichproben aus dieser Superpopulation.
    Dieses Model beruht aber auf der Annahme, dass der datengenerierende Prozess über die Zeit stabil bleibt. Das gilt offensichtlich nicht für Wahlen! Es stimmen nicht dieselben Menschen mit den selben Präferenzen ab. Wenn das der Fall wäre, wären selbst kleinste Unterschiede in den Wahlergebnissen unwahrscheinlich. Tatsächlich ändern Wähler aber ihre Meinung, weshalb auch große Unterschiede in Wahlergebnissen nicht selten sind.

Im Folgenden tue ich aber normalerweise so, als handele es sich nicht um Daten der Gesamtbevölkerung, sondern um eine Stichprobe, um zu zeigen, wie ich in dieser (deutlich typischeren) Situation verfahren würde.

b) Multiple Regression ¶

Eine einfache Regression zeigt den Zusammenhang zwischen zwei Variablen. Durch einen Plot mit Punkten wie die drei oben legt sie die Linie, die den kleinsten „durchschnittlichen“ vertikalen Abstand zu allen Punkten hat. Man kann nicht den tatsächlichen Durchschnitt der Abstände nehmen, da die Regression dann einfach immer eine horizontale Linie durch den Durchschnitt der abhängigen Variable (auf der y-Achse) legte, da der durchschnittliche Abstand zum Durchschnitt natürlich Null ist. Stattdessen wird mit dem Durchschnitt der quadrierten Abstände gearbeitet; auf Englisch die residual sum of squares (RSS).
    Eine multiple Regression macht prinzipiell dasselbe in höheren Dimensionen. Für zwei unabhängige Variablen legt es die zweidimensionale Ebene mit der geringsten RSS durch die Daten, für drei Variablen die dreidimensionale Hyperebene mit der geringsten RSS usw.

Der Output der unteren Code-Zelle zeigt nun einige Zahlen einer multiplen Regression, in der alle sozialen Variablen benutzt werden, um den Erststimmenanteil der SPD vorauszusagen. Relevant ist der mittlere Block und die Zahl ganz rechts oben:
    R-squared oder etwas klarer „variance explained“ misst, wie viel der Varianz in der abhängigen Variable (Erststimmenanteil) von der Regression beseitigt wird. In diesem Fall bleiben nur noch ca. 5% der ursprünglichen Varianz übrig! Solch hohen $R^2$ Werte sind in den Sozialwissenschaften extrem selten und indizieren eher, dass irgendetwas nicht stimmt.
    Der Block in der Mitte listet die Variablen dem Namen nach auf. Erst „const,“ was der konstante Faktor in der linearen Gleichung ist – was $a$ in $y=a+bx$ ist –, dann „num_bevölkerung“ usw. In der Spalte „coef“ wird der Wert des Koeffizienten vor der Variable angezeigt. Z. B. Hat die Arbeitslosenquote ca. den Koeffizienten 5,1. Das heißt, wenn wir alles andere gleich lassen und die Arbeitslosenquote um eine Einheit (hier einen Prozentpunkt) erhöhen, erhöht sich im Durchschnitt der Erststimmenanteil der SPD um 5,1 Einheiten.
    Die einzigen anderen wirklich wichtigen Zahlen sind in der vierten Spalte unter „$P>\vert t\vert $.“ Das sind die P-Werte der Koeffizienten, und wenn die Daten eine Stichprobe wären, gäben sie an, wie wahrscheinlich es ist, einen wenigstens so großen Koeffizienten per Zufall zu erhalten, obwohl es in der Realität keinen Zusammenhang gibt, der Koeffizient also Null ist. Z. B. hätte der Koeffizient der Variable „num_ausländer“ einen P-Wert von 0,425. Es gäbe also eine sehr hohe Wahrscheinlichkeit, dass der Koeffizient eigentlich Null ist; die Signifikanzgrenze wird standardmäßig bei 0.05 gezogen.

model = sm.OLS(daten[y], sm.add_constant(Stadtteile.drop(['Häfen', 'Stadt Bremen']))).fit()
model.summary()

c) Außreißer? ¶

Da es so wenig Stadtteile gibt, produziert Python eine Warnung, dass der Test auf Kurtose nicht funktioniert (der Computer kann nicht errechnen, ob die Fehler normalverteilt sind). Eine alternative ist sich einfach mit bloßem Auge den Plot der Residuen anzugucken (dem Abstand der Punkte von der Regressionslinie). Wenn sie nicht überall etwa den gleichen Abstand zu Null haben, sind die Daten heteroskedastisch oder es gibt einen nicht-linearen Zusammenhang zwischen den unabhängigen und der abhängigen Variable. Der untere Plot zeigt „studentized residuals,“ die gleichzeitig auch Ausreißer markieren. Ausreißer gibt es keine, allerdings könnte leichte Heteroskedastizität vorliegen; links gibt es schlicht zu wenig punkte, um das zu beurteilen.

studentized_residuals = model.get_influence().resid_studentized_internal

ax = subplots(figsize=(4, 4))[1]
ax.scatter(model.fittedvalues, studentized_residuals, s=5)

ax.scatter(model.fittedvalues[abs(studentized_residuals) > 3],
           studentized_residuals[abs(studentized_residuals) > 3],
           color='red', s=20, label='|Studentized Residual| > 3')

ax.set_xlabel('Fitted value')
ax.set_ylabel('Studentized Residual')
ax.axhline(0, c='k', ls='--')
ax.axhline(3, c='r')
ax.axhline(-3, c='r')
ax.legend(loc='upper left', bbox_to_anchor=(1, 1))

<matplotlib.legend.Legend at 0x13d636650>

d) Variance Inflation Factor ¶

Die Zusammenfassung der Regression oben sagt in der zweiten Fußnote auch, dass eventuell Multikollinearität vorliegt, also dass die unabhängigen Variablen stark untereinander korrelieren könnten.
    Es ist möglich, sich die paarweisen Korrelationen aller Variablen in einer Matrix anzeigen zu lassen, doch das ignoriert die Situation, dass Variable $X$ zwar nicht besonders stark mit den Variablen $A$, $B$ und $C$ korreliert, aber schon mit einer linearen Kombination der drei. Das bessere Maß ist der „variance inflation factor.“ Er misst wie stark jede einzelne Variable mit allen anderen Variablen gleichzeitig korreliert (statt nur paarweisen Korrelationen).
    Der Faktor kann nie unter 1 sein, ein Wert von 10 indiziert starke Multikollinearität, ein Wert von 5 stellt eine konservativere Schwelle dar … diese Nuancen sind aber völlig irrelevant für die Daten:

X = sm.add_constant(Stadtteile.drop(['Häfen', 'Stadt Bremen']))

vals = [VIF(X, i) for i in range(1, X.shape[1])]
vif = pd.DataFrame({'vif':vals}, index=X.columns[1:])
vif

Eine low-tech Herangehensweise wäre, einfach die Variablen mit dem höchsten VIF-Wert aus dem Modell zu nehmen, bis alle Werte unter 10 liegen. Da man den VIF des konstanten Faktors ignorieren kann, führt das zu diesem Modell:

X_reduced = X.drop(['num_ausländer', 'sgbII_pro_kopf', 'sgbII_leistungen_vor_wahl', 'anteil_erstes_quint'], axis=1)

vals = [VIF(X_reduced, i) for i in range(X_reduced.shape[1])]
vif = pd.DataFrame({'vif':vals}, index=X_reduced.columns)
vif

Wenn es sich um eine Stichprobe handelte, stünde aber noch nicht fest, dass alle diese Variablen auch statistisch signifikant sind. Das tue ich im nächsten Schritt.

e) ANOVA ¶

ANOVA steht für „analysis of variance.“ Eine Anwendung für ANOVA ist es, alle Variablen eines Modells einem Hypothesentest zu unterziehen: Um zu checken, ob die Variable „anteil_ausländer“ signifikant ist, vergleicht ANOVA den $R^2$-Wert des Modells mit der Variable und des Modells ohne sie. Wenn die Erhöhung des $R^2$-Wertes relativ zur Varianz der Residuen groß ist, ist die Variable signifikant.
    Einfacher: Wenn der P-Wert in der Spalte ganz rechts kleiner als 0.05 ist, ist die Variable wahrscheinlich nicht irrelevant.

# Renaming the first so that smf.ols does not complain.  
# Renaming the second for convenience later on.

daten_ = daten.copy()
daten_.rename(columns={'SPD - Erststimmenanteil (%)  2021': 'SPD_erststimmen_2021', 
                       'Wahlbeteiligung (%) 2021': 'wahlbeteiligung_2021'}, inplace=True)

formula = ' + '.join(X_reduced.columns)
model_test = smf.ols(f'SPD_erststimmen_2021 ~ {formula}', data=sm.add_constant(daten_)).fit()

anova_lm(model_test)

Es bleibt also nur der konstante Term, „anteil_ausländer“ und „anteil_letztes_quint;“ die letzte Variable bezieht sich auf das Fünftel mit dem höchsten Einkommen.
Nun kann ich also eine Regression mit diesem kleineren Modell errechnen:

# Fit the model without insignificant predictors
formula_reduced = 'SPD_erststimmen_2021 ~ anteil_ausländer + anteil_letztes_quint'
model_reduced = smf.ols(formula_reduced, data=daten_).fit()
print(model_reduced.summary())

                             OLS Regression Results                             
================================================================================
Dep. Variable:     SPD_erststimmen_2021   R-squared:                       0.686
Model:                              OLS   Adj. R-squared:                  0.644
Method:                   Least Squares   F-statistic:                     16.38
Date:                  Wed, 28 Aug 2024   Prob (F-statistic):           0.000169
Time:                          22:40:47   Log-Likelihood:                -46.606
No. Observations:                    18   AIC:                             99.21
Df Residuals:                        15   BIC:                             101.9
Df Model:                             2                                         
Covariance Type:              nonrobust                                         
========================================================================================
                           coef    std err          t      P>|t|      [0.025      0.975]
----------------------------------------------------------------------------------------
Intercept               46.8742      8.532      5.494      0.000      28.690      65.059
anteil_ausländer        -3.1939     23.710     -0.135      0.895     -53.730      47.342
anteil_letztes_quint   -59.8935     19.244     -3.112      0.007    -100.912     -18.875
==============================================================================
Omnibus:                        2.487   Durbin-Watson:                   1.612
Prob(Omnibus):                  0.288   Jarque-Bera (JB):                1.337
Skew:                          -0.666   Prob(JB):                        0.512
Kurtosis:                       3.099   Cond. No.                         38.4
==============================================================================

Notes:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.

Leider ist der P-Wert von „anteil_ausländer“ extrem hoch (0,895). Aber wenn ich nochmal ANOVA auf das kleinere Modell loslasse (Zelle direkt unten), sind beide Variablen signifikant.
    Das könnte daran liegen, dass es immer noch Multikollinearität im Modell gibt, da sie die Standardfehler erhöhte und die statistische „Power“ reduzierte. Aber wenn ich nochmal die variance inflation factors errechne (eine Zelle weiter unten), sind sie nicht besonders hoch … Diese low-tech Herangehensweise per hand signifikante Zusammenhänge herauszufinden, führt offensichtlich zu nichts. Ich muss andere Methoden benutzen.

print('ANOVA results:')
print(anova_lm(model_reduced)['PR(>F)'])

ANOVA results:
anteil_ausländer        0.000232
anteil_letztes_quint    0.007136
Residual                     NaN
Name: PR(>F), dtype: float64

# VIF without insignificant predictors
X_reduced_more = X.drop(['num_ausländer', 'sgbII_leistungen_vor_wahl', 'anteil_erstes_quint', 'arbeitslosenquote_vor_wahl', 'num_bevölkerung'], axis=1)

vals = [VIF(X_reduced, i) for i in range(X_reduced_more.shape[1])]
vif = pd.DataFrame({'vif':vals}, index=X_reduced_more.columns)
print(vif)

                             vif
const                 136.799558
anteil_ausländer        1.073771
anteil_letztes_quint    7.907788
sgbII_pro_kopf          4.104336

Faktorenanalyse ¶

Hoch

Da die einfachsten Methoden, die Daten zu analysieren, zu nur wenig einsichtsreichen Schlussfolgerungen geführt haben, benutze ich nun zwei moderat kompliziertere.
    Da Faktorenanalyse latente Faktoren aufspüren soll, die alle anderen Variablen beeinflusst, ergibt es Sinn, erstmal intuitiv darüber nachzudenken, was dieser Faktor bedeuten könnte. Die Variablen sind:

1. Anteil Ausländer.
2. Anteil im ärmsten Fünftel.
3. Arbeitslosenquote.
4. Wahlbeteiligung.
5. SGB-II Leistungen pro Kopf.
6. SGB-II Leistungen insgesamt.
7. Anzahl Ausländer insgesamt.
8. Anzahl Bevölkerung insgesamt.
9. Anteil im reichsten Fünftel.

Die drei „insgesamt“ Variablen sind v. a. Maße der Bevölkerungsanzahl. Es gibt keinen analytischen Grund anzunehmen, dass die Anzahl der Menschen in einem Stadtteil unabhängig von der Bevölkerungsdichte irgendeinen Einfluss auf das Wahlergebnis haben sollte. Wenn zwei absolut gleiche Stadtteile aufeinmal zusammengelegt würden, ginge die Bevölkerungszahl nach oben aber wir würden keine Veränderung im Wahlverhalten erwarten.
    Es bleiben also die anteilsmäßigen Variablen. Jedoch scheint mir intuitiv der Anteil im reichsten Fünftel qualitativ anders zu sein als die ersten fünf Variablen (es wird sich tatsächlich gleich herausstellen, dass diese fünf Variablen die beste Kombination sind).
    Was verbindet diese Variablen? Vielleicht so etwas wie „Sozialschwäche“ oder eine „prekäre Lage.“ Der genaue Name spielt für die Analyse keine Rolle. Aber es scheint wengstens erstmal so, als könnte es einen unterliegenden Faktor geben.

a) Sind die Daten für Faktorenanalyse geeignet? ¶

Auf bloße Spekulation und Eindrücke sollte man sich nicht verlassen. Als nächstes überprüfe ich daher, ob die Daten tatsächlich für Faktorenanalyse geeignet sind.

# Checking whether the data are even suitable for factor analysis

X_clean = pd.merge(X.drop(columns=['const', 
                                   'anteil_letztes_quint', 
                                   'num_ausländer', 
                                   'sgbII_leistungen_vor_wahl', 
                                   'num_bevölkerung']), 
                   daten_['wahlbeteiligung_2021'] , 
                   on='region')

kmo_all, kmo_model = calculate_kmo(X_clean)
print('KMO Test: ', round(kmo_model, 3))

bartlett_test = calculate_bartlett_sphericity(X_clean)
print('Bartlett’s test: ', bartlett_test[1])
print()
print('Benutzte Variablen:')
print(', '.join(X_clean.columns))

KMO Test:  0.872
Bartlett’s test:  1.4428985414960861e-24

Benutzte Variablen:
anteil_ausländer, anteil_erstes_quint, arbeitslosenquote_vor_wahl, sgbII_pro_kopf, wahlbeteiligung_2021

Das Kaiser-Meyer-Olkin Kriterium misst den Anteil der Varianz unter den Variablen, die allen Variablen gemeinsam sein könnte. Wenn also zwei Variablen miteinander aber nicht mit den anderen Variablen korrelieren, zählt das Kriterium das nicht. Wir stellen uns ja vor, dass es gemeinsame Faktoren gibt, die alle Variablen beeinflussen. Wenn also nicht alle Variablen mit allen anderen korrelieren, gibt es solche Faktoren eher nicht.
    KMO produziert Werte in $[0,1]$ und sollte nahe bei Eins sein; unter 0,8 gilt als mittelmäßig und unter 0,6 spricht gegen eine Faktorenanalyse.
    Bartletts Test produziert u. a. einen P-Wert, der nahe an Null sein sollte. Das Ergebnis des Tests hier hat 23 Nachkommastellen vor der ersten Eins. Die Notation e-24 ist dasselbe wie $\times 10^{-24}$.

Die Kombination aus Variablen ist letztendlich durch trial-and-error zustande gekommen geleitet von er Überlegung oben, wie der Faktor aussehen müsste. Einen höheren Werte für KMO gibt für keine andere Kombination.

b) Wie viele Faktoren? ¶

Eine etwas analytischere Methode, um zu entscheiden, wie viele Faktoren es gibt, ist es, den sogenannten „Scree“-Plot anzugucken:

# Performing Factor Analysis
fa = FactorAnalyzer()
fa.fit(X_clean)

ev, v = fa.get_eigenvalues()

# Scree Plot
plt.figure(figsize=(6, 5))
plt.plot(ev, 'o-')
plt.xlabel('Anzahl der Faktoren')
plt.ylabel('Eigenwert')
plt.title('Scree Plot (Faktorenanalyse)')

variables = '\n'.join([
    r'$\bf{Variablen\ (alles\ pro\ Stadtteil)\ }$',
    '',
    '1.) Ausländer als Anteil der Bevölkerung.',
    '2.) Anteil Bevölkerung im untersten Einkommensquintil.',
    '3.) Arbeitslosenquote im Jahr vor der Wahl.',
    '4.) SGB II Leistungen pro Kopf.',
    '5.) Wahlbeteiligung 2021.',
])

plt.figtext(0.5, -0.1, variables, ha="center", fontsize=10, bbox={"facecolor":"lightgrey", "alpha":0.5, "pad":5})
plt.subplots_adjust(bottom=0.25)

plt.show()

Man schaut, wo der „Ellenbogen“ im Graph ist, also der Punkt, nach dem der Plot flach wird.^† Die y-Achse mag etwas kryptisch wirken, sie misst prinzipiell, wie viel der Varianz erklärt wird. Anscheinend erklärt der erste Faktor fast die gesamte Varianz!
Wie viel er genau erklärt, wird in der nächsten Zelle ausgerechnet:

# Choose the number of factors based on the scree plot
num_factors = 1

fa_ = FactorAnalyzer(n_factors=num_factors)
fa_.fit(X_clean)

# Loadings and variance
factor_loadings = fa_.loadings_
print('Factor Loadings:\n', factor_loadings)

explained_variance = fa_.get_factor_variance()[1]
print("Explained Variance: ", explained_variance)

Factor Loadings:
 [[-0.95876235]
 [-0.94796932]
 [-0.9794939 ]
 [-0.98996176]
 [ 0.88025438]]
Explained Variance:  [0.90643029]

Die „Factor Loadings“ zeigen, wie stark die einzelnen Variablen mit dem ersten Faktor korrelieren. Die Wahlbeteiligung, die letzte Variable, hat als einzige eine positive Korrelation, was intuitiv Sinn ergibt: Mehr SGB-II Leistungen pro Kopf bedeuten eine prekärere Lage aber dasselbe gilt für weniger Wahlbeteiligung. Das Vorzeichen muss wechseln.
    Die „Explained Variance“ ist genau das, was der Name sagt. Anscheinend erklärt der erste Faktor 91% der Varianz in den fünf Variablen.

c) Regression mit Factor Scores ¶

Die fünf Variablen sind aber unsere unabhängigen Variablen und wir wollen ja die abhängige Variable, den Erststimmenanteil der SPD, erklären. Dafür sind die „Factor Scores“ nützlich. Sie zeigen für jeden Stadtteil an, wie stark der latente Faktor, also die Sozialschwäche oder Prekarität, zum Ausdruck kommt. Mit diesen Scores kann man dann wieder eine einfache Regression errechnen:

# Getting the factor scores
factor_scores = fa_.transform(X_clean)
X_factor = pd.DataFrame({'Factor Score':factor_scores.flatten()}, index=X_clean.index) 

# Fitting the regression model
X_factor = sm.add_constant(X_factor)
y_factor = daten_['SPD_erststimmen_2021']
model_factor = sm.OLS(y_factor, X_factor).fit()

print(model_factor.summary())

                             OLS Regression Results                             
================================================================================
Dep. Variable:     SPD_erststimmen_2021   R-squared:                       0.602
Model:                              OLS   Adj. R-squared:                  0.577
Method:                   Least Squares   F-statistic:                     24.20
Date:                  Sat, 24 Aug 2024   Prob (F-statistic):           0.000154
Time:                          18:00:40   Log-Likelihood:                -48.738
No. Observations:                    18   AIC:                             101.5
Df Residuals:                        16   BIC:                             103.3
Df Model:                             1                                         
Covariance Type:              nonrobust                                         
================================================================================
                   coef    std err          t      P>|t|      [0.025      0.975]
--------------------------------------------------------------------------------
const           32.7353      0.907     36.090      0.000      30.812      34.658
Factor Score    -4.4871      0.912     -4.919      0.000      -6.421      -2.553
==============================================================================
Omnibus:                        2.571   Durbin-Watson:                   1.679
Prob(Omnibus):                  0.277   Jarque-Bera (JB):                1.045
Skew:                          -0.540   Prob(JB):                        0.593
Kurtosis:                       3.475   Cond. No.                         1.01
==============================================================================

Notes:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.

Der $R^2$-Wert ganz rechts oben zeigt an, dass der Faktor etwa 60% der Varianz in den Erststimmenanteilen der SPD zwischen den Stadtteilen erklärt. Die P-Werte sind bedeutungslos, aber wenn es sich um eine Stichprobe handelte, wäre der Effekt des Faktors hochsignifikant, da p < 0.0005.

d) Signifikanz ¶

„Signifikanz“ bedeutet in diesem und allen anderen Zusammenhängen nie, dass ein Koeffizient oder ein Unterschied praktisch relevant oder bedeutend ist. Signifikanz bezieht sich stets nur auf die Wahrscheinlichkeit, dass der Wert eines Koeffizienten per Zufall Zustande gekommen wäre, wenn es eigentlich keinen Zusammenhang gibt. „Hochsignifikant“ bedeutet daher, dass es eine sehr kleine Wahrscheinlichkeit gibt, einen so extremen Wert zu bekommen, obwohl es eigentlich keinen Zusammenhang gibt.
    Wenn es in einer Stichprobe einen Unterschied im Notendurchschnitt der Abiture von Männern und Frauen von 0,004 Punkten gibt – z. B. 2,1 und 2,104 –, könnte dieser völlig irrelevante Unterschied trotzdem hochsignifikant sein, wenn die Stichprobe groß genug ist. Wenn es diesen Unterschied im Durchschnitt aller Abiture statt nur in einer Stichprobe gibt, ist das Konzept von „Signifikanz“ nicht anwendbar, da keine Unsicherheit in unseren Daten ist.

e) Auswertung ¶

Auch wenn das Modell nur auf den Erststimmen trainiert wurde, ist seine Performance für den Zweitstimmenanteil nicht viel schlechter, wie der Plot unten zeigt. In welche Richtung die Kausalität geht, machen die Daten aber nicht klar:
    Ist es, dass politisches Verhalten von SPD-Wählern v. a. von diesem Faktor beeinflusst wird und es daher kaum einen Unterschied in Erst- und Zweitstimmen gibt? Oder gibt es aus unabhängigen Gründen kein Stimmensplitting, weshalb der Faktor zufällig beide Stimmenanteile erklärt?

Wie dem auch sei, zeigen die Tabellen unter dem Plot, wie viel der Faktor Sozialschwäche (oder prekäre Lage oder wie man ihn nennen möchte) von den Stimmenunterschieden zwischen den Stadtteilen für verschiedene Parteien erklärt. Die relevante Zeile ist jeweils die letzte, die „R^2“ heißt.
    Etwa 60% werden für SPD- und AfD-Erstimmen erklärt und etwa ein Viertel für Grüne und FDP aber mit anderem Vorzeichen! Das Vorzeichen kann man dem $R^2$ natürlich nicht ansehen, aber zwei Zeilen drüber ist die Korrelation, die auch negativ sein kann. Je sozial schwächer ein Stadtteil, desto höher der SPD und AfD Stimmenanteil aber desto niedriger der FDP und Grünen Stimmenanteil.
    Das bedeutet jedoch nicht, dass sich SPD- und AfD-Wähler ähnlich wären oder für einander substituiert werden könnten. Was analysiert wird sind die Stadtteile, in denen sie leben. Es könnte also sein, dass es (wenigstens) zwei Typen von Wählern gibt, die jeweils unterschiedlich auf höhere Prekarität reagieren oder sie verursachen; SPD- und AfD-Wähler könnten sich beide von sozialer Schwäche angezogen fühlen oder sie könnten beide auf ähnliche Art sozial inkompetent sein und daher die Prekarität ihrer Stadtteile erzeugen.
    Dass es bei der Linken keinen starken Zusammenhang zwischen dem Faktor und ihren Stimmen gibt, könnte auch an zwei unterliegenden, entgegengesetzten Tendenzen liegen: Leute in sehr prekären und in sehr unprekären Stadtteilen wählen sie, während der Effekt für die anderen Parteien homogener ist.

Diese Korrelationen zeigen auch nicht an, mit welcher Partei der Wettkampf um Stimmen besonders hoch ist. Wenn SPD- und AfD-Wähler tatsächlich nicht für einander substituiert werden können, bedeutet eine Stimme mehr für die AfD nicht direkt eine weniger für die SPD. Und auch wenn Grünen- und SPD-Wähler anders auf soziale Schwäche reagieren, könnte es trotzdem in jedem einzelnen Stadtteil sein, dass eine Stimme mehr für die Grünen eine weniger für die SPD bedeutete.
Dieser Frage widme ich mich im nächsten Teil.

plt.figure(figsize=(6, 6))
plt.scatter(X_factor['Factor Score'], y_factor, label='Erstimmenanteil', color='blue')
plt.scatter(X_factor['Factor Score'], daten['SPD - Zweitstimmenanteil (%)  2021'], marker='^', label='Zweitstimmenanteil', color='darkorange')

plt.xlabel('Prekäre Lage (negativ = prekärer)')
plt.ylabel('SPD Stimmenanteil')

plt.legend(loc='best')


regions = X_factor.index
for i, region in enumerate(regions):
    plt.text(X_factor['Factor Score'][i], y_factor[i], region, fontsize=8, ha='right')

plt.show()

parties = ['SPD', 'CDU', 'GRÜNE', 'DIE LINKE', 'FDP', 'AfD']
result_types = ['Coefficient', 'Correlation', 'P-Value', 'R^2']

Parteien_factor_erststimmen = pd.DataFrame(index=result_types)
Parteien_factor_zweitstimmen = pd.DataFrame(index=result_types)

for ordinal in ['Erst', 'Zweit']:
    for party in parties:
        y_party = daten[f'{party} - {ordinal}stimmenanteil (%)  2021']
        model_party = sm.OLS(y_party, X_factor).fit()
        
        coefficient_party = model_party.params['Factor Score']
        p_value_party = model_party.pvalues['Factor Score']
        correlation_party, _ = pearsonr(X_factor['Factor Score'], y_party)
        rsquared = correlation_party ** 2

        if ordinal == 'Erst':
            Parteien_factor_erststimmen[party] = [
                coefficient_party,
                correlation_party,
                p_value_party,
                rsquared]
        else:
            Parteien_factor_zweitstimmen[party] = [
                coefficient_party,
                correlation_party,
                p_value_party,
                rsquared]

Parteien_factor_erststimmen = Parteien_factor_erststimmen.round(3)
Parteien_factor_zweitstimmen = Parteien_factor_zweitstimmen.round(3)

Parteien_factor_erststimmen

Parteien_factor_zweitstimmen

LASSO und Validierung ¶

Hoch

Formal ist der unterschied zwischen LASSO und normaler Regression, dass die Koeffizienten in einer LASSO-Regression folgenden Term minimieren: $$ \text{RSS} + \lambda\sum^p_{j=1}\lvert \beta_j \rvert $$ Die $\beta_j$ sind die Koeffizienten, von denen es so viele gibt, wie unsere Tabelle Spalten hat. Eine „Tabelle“ in der Statistik heißt eigentlich Designmatrix, und wenn man die im Allgemeinen behandelt, sagt man normalerweise es sei eine $n\times p$ Matrix, also eine Matrix mit $n$ Zeilen und $p$ Spalten. Daher kommt das $p$ auf dem Summenzeichen. „RSS“ ist immer noch die „residual sum of squares,“ die den Abstand aller Punkte von der Regressionslinie misst.
    Eine normale Regression minimiert einfach die RSS. Aber LASSO fügt dem noch die Gesamtgröße der Koeffizienten hinzu; je mehr Koeffizienten es gibt und je größer sie sind, desto kleiner müsste die RSS sein, um das auszugleichen. Der Faktor $\lambda$ kontrolliert, wie stark das Modell für viele und große Koeffizienten bestraft wird (im Englischen tatsächlich „penalty term“).
    Ein letztes Detail ist, dass der absolute Betrag der Koeffizienten in die Summe kommt. Man könnte auch die quadrierten Koeffizienten nehmen, dann wäre es eine Ridge-Regression, aber mit dem absoluten Betrag bekommt LASSO einige Koeffizienten auf genau Null, statt nur sehr klein. Deshalb ist LASSO in der Lage, die wichtigsten variablen auszusuchen, während z. B. Ridge fast immer alle Variablen im Modell lässt. LASSO funktioniert daher am besten, wenn tatsächlich nur eine Handvoll Variablen wichtig sind.

a) Das LASSO-Modell anpassen ¶

Als nächstes passe ich eine LASSO-Regression auf die gesamten Daten an: Die sozialen Daten und auch die Stimmenergebnisse der verschiedenen Parteien und deren Veränderung von 2017 auf 2021.

Zwei technische Besonderheiten; kann ohne Verständnisverlust übersprungen werden:     Da für LASSO die Größe der Koeffizienten eine Rolle spielt, sollten die Daten normalisiert werden. Wenn ich militärische Macht mit einer Regression messen wollte und Anzahl der Soldaten und Militärausgaben in Euro als Variablen hätte, würde die erste Variable Zahlen von mehreren Hunderttausend annehmen aber die zweite von mehreren Milliarden. Normaler Regression ist das egal aber LASSO könnte die Ausgaben unterbewerten. Nach einer Normalisierung werden alle Variablen in Standardabweichungen von ihrem Mittelwert gemessen.     Zweitens ist das $\lambda$ ein Hyperparameter, der getunt werden muss. LASSO weiß nicht automatisch, was der beste $\lambda$-Wert ist. Das passiert hier dadurch, dass ich den `LassoCV`-*estimator* aus `scikit-learn` benutze statt des normalen `Lasso`. „CV“ steht für „*cross validation*,“ eine äußerst clevere und simple Methode, um ohne einen extra Datensatz zum Testen des Modells auszukommen: Die Daten werden in z. B. Fünftel geteilt, und dann werden immer vier Teile benutzt, um das Modell zu trainieren, und der fünfte, um das Modell zu testen. Das macht man mit allen fünf Teilen und nimmt den Durchschnitt der Ergebnisse.

daten.columns
daten_lasso = daten.drop(columns=['SPD - Erststimmenanteil (%)  2021', 
                                  'SPD - Erststimmenanteil: Veränderung 2021-2017', 
                                  'SPD - Zweitstimmenanteil (%)  2021', 
                                  'SPD - Zweitstimmenanteil: Veränderung 2021-2017'])
y_lasso = daten['SPD - Erststimmenanteil (%)  2021']

# Pipeline for LASSO
lasso_cv = LassoCV(cv=5, max_iter=10000)
pipeline = Pipeline([
    ('scaler', StandardScaler()),
    ('lasso_cv', lasso_cv)
])

# Fit the model
pipeline.fit(daten_lasso, y_lasso)

# Get the best alpha (penalty term) and the corresponding model
best_alpha = pipeline.named_steps['lasso_cv'].alpha_
best_model = pipeline.named_steps['lasso_cv']

# Uncomment as necessary. Cluttering up my analysis right now.
'''
print(f"Best alpha (penalty term): {best_alpha}")
print("Coefficients:", best_model.coef_)
'''

# If you want to evaluate the performance:
y_pred = pipeline.predict(daten_lasso)
rmse = np.sqrt(mean_squared_error(y_lasso, y_pred))
print(f"Root Mean Squared Error: {rmse:.3f}")

Root Mean Squared Error: 0.777

Der „root mean square error“ direkt oben ist fast das, was der Name sagt: „root“ und „square“ heben sich hier im Prinzip auf und die 0,777 gibt an, dass das LASSO-Modell im Schnitt soviele Prozentpunkte vom tatsächlichen Wert der Erststimmenanteile in den Stadtteilen entfernt ist. Das Modell kommt also sehr nahe an die Realität heran!
    Nun ist es Zeit, die tatsächlichen Koeffizienten anzusehen. Sie haben im Gegensatz zur normalen Regression keine intrinsische Bedeutung, da alle Variablen normalisiert wurden. Aber ihre relative Größe indiziert die Wichtigkeit der Variable und ihr Vorzeichen sagt auch, in welche Richtung die Variable mehr Stimmen bedeutet:

coef_df = pd.DataFrame({
    'Variable': daten_lasso.columns,
    'Coefficient': best_model.coef_
})

#print(lasso_cv.intercept_)
print(coef_df)

                                             Variable  Coefficient
0                                     num_bevölkerung     0.264051
1                                       num_ausländer     0.000000
2                                    anteil_ausländer     0.076958
3                                 anteil_erstes_quint     0.491252
4                                anteil_letztes_quint    -0.849459
5                          arbeitslosenquote_vor_wahl     0.000000
6                           sgbII_leistungen_vor_wahl    -0.000000
7                                      sgbII_pro_kopf    -0.000000
8                            Wahlbeteiligung (%) 2021    -0.000000
9              Wahlbeteiligung: Veränderung 2021-2017    -0.983569
10                  CDU - Erststimmenanteil (%)  2021    -0.000000
11     CDU - Erststimmenanteil: Veränderung 2021-2017    -0.000000
12            DIE LINKE - Erststimmenanteil (%)  2021    -1.200486
13  DIE LINKE - Erststimmenanteil: Veränderung 202...     0.776733
14                GRÜNE - Erststimmenanteil (%)  2021    -0.000000
15   GRÜNE - Erststimmenanteil: Veränderung 2021-2017    -0.000000
16                  AfD - Erststimmenanteil (%)  2021     0.000000
17     AfD - Erststimmenanteil: Veränderung 2021-2017    -0.000000
18                  FDP - Erststimmenanteil (%)  2021    -0.000000
19     FDP - Erststimmenanteil: Veränderung 2021-2017    -0.000000
20                 CDU - Zweitstimmenanteil (%)  2021    -0.000000
21    CDU - Zweitstimmenanteil: Veränderung 2021-2017    -0.000000
22           DIE LINKE - Zweitstimmenanteil (%)  2021    -0.000000
23  DIE LINKE - Zweitstimmenanteil: Veränderung 20...     0.000000
24               GRÜNE - Zweitstimmenanteil (%)  2021    -1.316454
25  GRÜNE - Zweitstimmenanteil: Veränderung 2021-2017    -0.000000
26                 AfD - Zweitstimmenanteil (%)  2021     0.000000
27    AfD - Zweitstimmenanteil: Veränderung 2021-2017    -1.080452
28                 FDP - Zweitstimmenanteil (%)  2021    -2.491543
29    FDP - Zweitstimmenanteil: Veränderung 2021-2017     0.000000

b) Interpretation des LASSO-Modells ¶

Es wird sich gleich herausstellen, dass diese Koeffizienten leider nicht besonders stabil sind. Aber exemplum gratia tue ich so, als seien sie stabil, und erkläre, wie man diese Koeffizienten interpretieren könnte:

1.) Die Wahlbeteiligung hat einen Nullkoeffizienten, aber die Veränderung der Wahlbeteiligung hat einen moderat negativen.
    Diese Kombination indiziert, dass der Wählerblock der SPD stabiler ist, als der anderer Parteien. Wenn sich die Wahlbeteiligung erhöht und die rohe Anzahl der SPD-Wähler nach oben geht, muss sich die Anzahl anderer Parteien noch mehr erhöhen, da ein Anstieg in der Wahlbeteiligung ja zu einem geringeren Anteil für die SPD führt. Das könnte daran liegen, dass die Wählerschaft der anderen Parteien wankelmütiger oder leichter zu mobilisieren ist.

2.) Der Erststimmenanteil der Linken hat einen solide negativen Koeffizienten, aber die Veränderung dieses Anteils einen moderat positiven.
    Was auch immer dazu führt, dass Wähler der Linken zur Wahl gehen, führt auch dazu, dass SPD-Wähler zur Wahl gehen (der positive Koeffizient). Aber unter sonst gleichen Bedingungen, einschließlich der Gesamtwahlbeteiligung (!), bedeutet eine Stimme mehr für die Linke eine Stimme weniger für die SPD, da es direkte Konkurrenz gibt (der negative Koeffizient). Es könnte auch sein, dass die Linke, unter sonst gleichen Bedingungen, Stimmen von der SPD stiehlt, aber noch mehr von der Konkurrenz stiehlt; dann könnte ein erhöhter Stimmenanteil für die Linke weniger für die SPD bedeuten, aber „noch mehr weniger“ für die Konkurrenz außerhalb von Linke und SPD.

3.) Die Kombination der AfD-Zweitstimmen und deren Veränderung zeigt an, dass es keine direkte Konkurrenz um Wählerstimmen zwischen SPD und AfD gibt (Koeffizient für Stimmenanteil ist Null). Und AfD-Wähler werden durch andere – oder eher gegensätzliche – Dinge mobilisiert als SPD-Wähler, die weniger erscheinen, wenn die AfD mehr mobilisiert – der negative Koeffizient vor der Veränderung der Zweitstimmen –; das bedeutet nicht unbedingt, dass absolut weniger SPD-Wähler erscheinen, sondern nur relativ weniger. Was auch immer einen größeren Anstieg der Wahlbeteiligung der AfD verursacht (also ein Anstieg von 0 auf 30 statt von 20 auf 30), führt auch zu einer geringeren Wahlbeteiligung für die SPD oder zu mehr Stimmen für andere Parteien.

4.) Die Zweitstimmenanteile für FDP und Grüne haben beide stark negative Koeffizienten. Die Veränderung dieser Anteile hat jeweils Nullkoeffizienten.
    Es gibt also einen starken Wettbewerb zwischen den Grünen/FDP und der SPD, was bedeutet, dass in Stadtteilen, in denen die Grünen oder die FDP einen höheren Stimmenanteil haben, die SPD tendenziell einen niedrigeren Anteil hat. Die Faktoren, die den Stimmenanteil der Grünen oder der FDP im Laufe der Zeit erhöhen oder verringern, haben jedoch keinen bedeutenden Einfluss auf den Stimmenanteil der SPD. (Das ist das intuitive Ergebnis. Warum sollte die Veränderung statt der tatsächliche Stimmenanteil irgendeine Rolle spielen?)

Jetzt passe ich ein LASSO-Modell nur auf die Stimmenergebnisse an:

daten_lasso_stimmen = daten_lasso.iloc[:, 8:]

pipeline.fit(daten_lasso_stimmen, y_lasso)
y_pred_stimmen = pipeline.predict(daten_lasso_stimmen)
rmse_stimmen = np.sqrt(mean_squared_error(y_lasso, y_pred_stimmen))
print(f"Root Mean Squared Error: {rmse_stimmen:.3f}")
best_model = pipeline.named_steps['lasso_cv']

coef_df_stimmen = pd.DataFrame({
    'Variable': daten_lasso_stimmen.columns,
    'Coefficient': best_model.coef_
})

print(coef_df_stimmen)

Root Mean Squared Error: 0.164
                                             Variable  Coefficient
0                            Wahlbeteiligung (%) 2021    -0.624041
1              Wahlbeteiligung: Veränderung 2021-2017     0.000000
2                   CDU - Erststimmenanteil (%)  2021    -1.520083
3      CDU - Erststimmenanteil: Veränderung 2021-2017    -0.104563
4             DIE LINKE - Erststimmenanteil (%)  2021    -2.522169
5   DIE LINKE - Erststimmenanteil: Veränderung 202...     0.557489
6                 GRÜNE - Erststimmenanteil (%)  2021    -0.000000
7    GRÜNE - Erststimmenanteil: Veränderung 2021-2017    -0.000000
8                   AfD - Erststimmenanteil (%)  2021    -0.000000
9      AfD - Erststimmenanteil: Veränderung 2021-2017    -0.000000
10                  FDP - Erststimmenanteil (%)  2021    -0.000000
11     FDP - Erststimmenanteil: Veränderung 2021-2017    -0.000000
12                 CDU - Zweitstimmenanteil (%)  2021    -2.239762
13    CDU - Zweitstimmenanteil: Veränderung 2021-2017     0.000000
14           DIE LINKE - Zweitstimmenanteil (%)  2021    -0.652327
15  DIE LINKE - Zweitstimmenanteil: Veränderung 20...    -0.581354
16               GRÜNE - Zweitstimmenanteil (%)  2021    -6.242644
17  GRÜNE - Zweitstimmenanteil: Veränderung 2021-2017    -0.418176
18                 AfD - Zweitstimmenanteil (%)  2021    -2.902884
19    AfD - Zweitstimmenanteil: Veränderung 2021-2017    -0.492824
20                 FDP - Zweitstimmenanteil (%)  2021    -2.667513
21    FDP - Zweitstimmenanteil: Veränderung 2021-2017    -0.342946

Der root mean square error ist mit 0,164 etwas besser als beim anderen Modell. Das wird daran liegen, dass in diesem Modell mehr Koeffizienten ungleich Null sind; das ist nichts gutes, da das Modell eher „überangepasst“ ist, wie sich gleich noch zeigt.
    Im allgemeinen sind die Koeffizienten nicht besonders robust. Viele sind größer/kleiner als beim ersten Modell oder haben sogar andere Vorzeichen. Man sollte daher nicht mit der Lupe über die Ergebnisse gehen (auch meine Interpretationen für die Koeffizienten oben sind daher in der Realität bedeutungslos). Grobkörnig zeigt sich etwa, dass bei den Zweitstimmen härtere Konkurrenz besteht – bis auf bei der Linken –, da bei den Erststimmen häufiger Koeffizienten gleich Null sind. Außerdem zeigt sich, dass die Anteile wichtiger sind als die Veränderungen (warum auch nicht? Genau so hätte man es erwartet). Außerdem gibt es keine großen Unterschiede in der Intensität der Konkurrenz zwischen den unterschiedlichen Parteien und der SPD, da die Koeffizienten immer etwa bei $-2,5$ liegen, bis auf den Ausreißer beim Zweitstimmenanteil der Grünen.

Die Wahlbeteiligung ist wohl relativ nicht so entscheidend, auch wenn die einfache Korrelation stark negativ ist (siehe ganz am Anfang der Analyse). Das impliziert, dass die einfache Korrelation durch Multikollinearität verursacht wird und nicht „echt“ ist. LASSO ist robuster gegenüber Multikollinearität und spürt korrekt auf, dass die stark negative Korrelation nicht unabhängig Informationen vermittelt.

c) Validierung der Modelle ¶

Der beste Test eines statistischen Modells ist, es auf einen frischen Datensatz loszulassen; auch der Grund, warum Replikation von Experimenten so wichtig ist, denn jedes neue Experiment produziert vorher noch nie dagewesene Daten. Ich könnte die Modelle an den Wahldaten 2023 (oder 2017) testen. Stattdessen wähle ich die typischere (weil einfachere) Variante, meinen Datensatz in zwei Teile aufzusplitten. Ein Fünftel wird mein Testdatensatz, auf dem rest werden die Modelle trainiert.
    Die zwei Modelle sind das volle Modell mit allen sozialen Daten und den Wahlergebnissen für die Parteien und das kleinere Modell mit nur den Wahlergebnissen:

X_train, X_test, y_train, y_test = train_test_split(daten_lasso, y_lasso, test_size=0.2, random_state=42)
X_train_, X_test_, y_train_, y_test_ = train_test_split(daten_lasso_stimmen, y_lasso, test_size=0.2, random_state=42)

# How well does the full model perform
pipeline.fit(X_train, y_train)

y_pred = pipeline.predict(X_test)
rmse_full = np.sqrt(mean_squared_error(y_test, y_pred))
print(f"Root Mean Squared Error, full model: {rmse_full:.3f}")

# How well does the smaller model perform
pipeline.fit(X_train_, y_train_)

y_pred_ = pipeline.predict(X_test_)
rmse_small = np.sqrt(mean_squared_error(y_test_, y_pred_))
print(f"Root Mean Squared Error, small model: {rmse_small:.3f}")

Root Mean Squared Error, full model: 1.751
Root Mean Squared Error, small model: 1.821

Obwohl oben das kleinere Modell einen kleineren Fehler hatte, als das volle Modell, hat es hier einen größeren, wenn ich es auf ungesehenen Daten teste. Das indiziert, dass das kleinere Modell überangepasst ist!
    In allen Daten ist immer auch etwas „noise,“ also zufällige Variation. Wenn statistische Modelle zu nah an ihre Trainingsdaten angepasst werden, tun sie quasi so, als sei dieses noise eigentlich Signal/Information. Wenn sie dann neue Daten voraussagen sollen, benutzen sie dafür z. T. noise, was zwangsläufig zu schlechteren Voraussagen führt (da noise eben per Definition nicht informativ ist).
    Wie kann es sein, dass ein kleineres Modell schlimmer überangepasst wird, als ein größeres? Wie oben erwähnt, hat das volle Modell tatsächlich weniger Koeffizienten, die ungleich Null sind. Daher ist das volle Modell effektiv kleiner als das „kleine“ Modell.

d) Verstehen vs. Voraussagen
An dieser Stelle lohnt es sich auch, auf den Unterschied zwischen Inferenz und Voraussagen aufmerksam zu machen. Das lineare Modell, das die multiple Regression am Anfang repräsentierte, war letztendlich nicht zu interpretieren, da die unabhängigen Variablen zu stark miteinander korrelierten. Es wurde nicht klar, was signifikant ist, oder welche Variablen man aus dem Modell raus lassen sollte. Das LASSO-Modell war deutlich einfacher zu interpretieren, da die Koeffizienten mehr oder weniger das bedeuteten, was sie zu sagen schienen.

Das bedeutet aber nicht, dass LASSO die besseren Voraussagen macht! In der nächsten Zelle teste ich die Voraussagen einer multiplen Regression, die genauso auf den vier Fünfteln der Daten trainiert wurde, auf dem letzten Fünftel:

LinReg = LinearRegression()
linreg_pred = LinReg.fit(X_train, y_train).predict(X_test)

rmse_linreg = np.sqrt(mean_squared_error(y_test, linreg_pred))
print(f"Root Mean Squared Error, linear regression: {rmse_linreg:.3f}")

Root Mean Squared Error, linear regression: 1.071

Diese Spannung zwischen Verstehen und Voraussagen (oder inference und prediction) besteht noch stärker für die mächtigsten Modelle im machine learning. Extreme Gradient Boosting oder Neural Networks sind völlig undurchsichtig, sind aber auch auf Daten, für die sie gut performen, unschlagbar! Wenn man ein boosting Modell angepasst hat, kann man sich keine Parameter ausdrucken lassen, die einem Einsicht in die unterliegende kausale Struktur geben. Aber da das Modell so sensibel gegenüber nicht-linearen Effekten und Interaktionen ist, kommt kein lineares Modell dagegen an, wenn die Realität denn tatsächlich aus vielen und komplizierten Nichtlinearitäten und Interaktionen besteht.

Mit Kanonen auf Spatzen: Random Forests und Boosting ¶

Hoch

Die Kurzzusammenfassung dieses Teils ist, dass stark nicht-lineare Modelle schlechter performen, als die linearen Modelle. Das indiziert, dass der Zusammenhang zwischen den Stimmenanteilen der SPD und den sozialen Daten und Wahlergebnissen der anderen Parteien grundsätzlich linear und unkompliziert ist.

Die Performance dieser Modelle könnte verbessert werden durch intensiveres Hyperparameter-Tuning. Im Gegensatz zur linearen Regression haben diese Modelle viele Parameter, die vom Charakter wie das $\lambda$ bei LASSO sind: Die Anzahl der „trees;“ die Tiefe der trees oder die Größe der Blätter; die Anzahl der Variablen, die in jedem Split benutzt wird; pruning oder kein pruning, wenn ja, wie stark; … All das sind variablen die durch cross validation optimiert werden müssen und nicht mathematisch a priori bestimmt werden können.

daten_tree = daten_lasso.copy()
y_tree = daten['SPD - Erststimmenanteil (%)  2021']

rf_model = RandomForestRegressor(random_state=42)

# RMSE for the RF regressor
rf_model.fit(X_train, y_train)
y_pred = rf_model.predict(X_test)

rmse_test = np.sqrt(mean_squared_error(y_test, y_pred))

# Feature importance
feature_importances = pd.DataFrame({'Feature': daten_tree.columns, 'Importance': rf_model.feature_importances_})
feature_importances = feature_importances.sort_values(by='Importance', ascending=False)

print(feature_importances)

                                              Feature  Importance
28                 FDP - Zweitstimmenanteil (%)  2021    0.085837
18                  FDP - Erststimmenanteil (%)  2021    0.084939
8                            Wahlbeteiligung (%) 2021    0.080414
16                  AfD - Erststimmenanteil (%)  2021    0.073383
17     AfD - Erststimmenanteil: Veränderung 2021-2017    0.071408
4                                anteil_letztes_quint    0.052340
2                                    anteil_ausländer    0.050784
9              Wahlbeteiligung: Veränderung 2021-2017    0.048676
6                           sgbII_leistungen_vor_wahl    0.044080
15   GRÜNE - Erststimmenanteil: Veränderung 2021-2017    0.042215
5                          arbeitslosenquote_vor_wahl    0.038396
27    AfD - Zweitstimmenanteil: Veränderung 2021-2017    0.031970
7                                      sgbII_pro_kopf    0.028558
19     FDP - Erststimmenanteil: Veränderung 2021-2017    0.027210
26                 AfD - Zweitstimmenanteil (%)  2021    0.026963
14                GRÜNE - Erststimmenanteil (%)  2021    0.025049
10                  CDU - Erststimmenanteil (%)  2021    0.024914
29    FDP - Zweitstimmenanteil: Veränderung 2021-2017    0.024019
25  GRÜNE - Zweitstimmenanteil: Veränderung 2021-2017    0.023754
3                                 anteil_erstes_quint    0.021371
21    CDU - Zweitstimmenanteil: Veränderung 2021-2017    0.020593
13  DIE LINKE - Erststimmenanteil: Veränderung 202...    0.011849
11     CDU - Erststimmenanteil: Veränderung 2021-2017    0.010453
22           DIE LINKE - Zweitstimmenanteil (%)  2021    0.010108
0                                     num_bevölkerung    0.009895
1                                       num_ausländer    0.009288
12            DIE LINKE - Erststimmenanteil (%)  2021    0.007886
20                 CDU - Zweitstimmenanteil (%)  2021    0.005101
23  DIE LINKE - Zweitstimmenanteil: Veränderung 20...    0.004771
24               GRÜNE - Zweitstimmenanteil (%)  2021    0.003777

Es gibt zwar keine Parameter, die man sich angucken könnte, aber es gibt „feature importances,“ die allerdings eine etwas abstrakte Größe messen. Wie man direkt oben sieht, gibt es auch kein richtiges Muster. Und die feature importances beim boosting (zwei Zellen weiter unten) sind relativ verschieden. Ich würde diese zahlen nicht zu ernst nehmen.

Die nächste Zelle zeigt, dass der random forest schlechter performt, als die linearen Modelle:

rf_model.fit(daten_tree, y_tree)
full_pred_rf = rf_model.predict(daten_tree)

rmse_notest = np.sqrt(mean_squared_error(y_tree, full_pred_rf))

print(f'Root Mean Squared Error, no validation: {rmse_notest:.4f}')
print(f'Root Mean Squared Error, yes validation: {rmse_test:.4f}')

Root Mean Squared Error, no validation: 1.1972
Root Mean Squared Error, yes validation: 4.0213

Der Fehler des Modells ist mehrmals größer als der Fehler vom LASSO, sowohl wenn ich den random forest auf frischen Daten teste, als auch wenn nicht.

Als letztes passe ich das stärkste mir bekannte Modell an: XGBoost. Man vergleiche die feature importances mit denen von oben und erkenne kein Muster. Darunter ist der Fehler des Modells als root mean square error.

xgb_reg = xgb.XGBRegressor(objective='reg:squarederror', random_state=42)
xgb_reg.fit(X_train, y_train)

# How well does the model perform
y_pred = xgb_reg.predict(X_test)
rmse_validation = mean_squared_error(y_test, y_pred, squared=False)

# Feature importances 
importances_boost = xgb_reg.feature_importances_
importance_df = pd.DataFrame({'Feature': X_train.columns, 'Importance': importances_boost})
print(importance_df.sort_values(by='Importance', ascending=False))

                                              Feature    Importance
16                  AfD - Erststimmenanteil (%)  2021  6.410742e-01
2                                    anteil_ausländer  1.543060e-01
18                  FDP - Erststimmenanteil (%)  2021  8.559006e-02
28                 FDP - Zweitstimmenanteil (%)  2021  8.475830e-02
9              Wahlbeteiligung: Veränderung 2021-2017  2.053055e-02
0                                     num_bevölkerung  5.657035e-03
1                                       num_ausländer  5.275896e-03
13  DIE LINKE - Erststimmenanteil: Veränderung 202...  2.740611e-03
14                GRÜNE - Erststimmenanteil (%)  2021  3.839945e-05
3                                 anteil_erstes_quint  1.297326e-05
5                          arbeitslosenquote_vor_wahl  7.901166e-06
6                           sgbII_leistungen_vor_wahl  6.134392e-06
4                                anteil_letztes_quint  1.325273e-06
21    CDU - Zweitstimmenanteil: Veränderung 2021-2017  4.358324e-07
11     CDU - Erststimmenanteil: Veränderung 2021-2017  3.147557e-08
19     FDP - Erststimmenanteil: Veränderung 2021-2017  2.206468e-08
8                            Wahlbeteiligung (%) 2021  2.025824e-08
29    FDP - Zweitstimmenanteil: Veränderung 2021-2017  1.504093e-08
17     AfD - Erststimmenanteil: Veränderung 2021-2017  0.000000e+00
12            DIE LINKE - Erststimmenanteil (%)  2021  0.000000e+00
20                 CDU - Zweitstimmenanteil (%)  2021  0.000000e+00
10                  CDU - Erststimmenanteil (%)  2021  0.000000e+00
22           DIE LINKE - Zweitstimmenanteil (%)  2021  0.000000e+00
23  DIE LINKE - Zweitstimmenanteil: Veränderung 20...  0.000000e+00
24               GRÜNE - Zweitstimmenanteil (%)  2021  0.000000e+00
25  GRÜNE - Zweitstimmenanteil: Veränderung 2021-2017  0.000000e+00
26                 AfD - Zweitstimmenanteil (%)  2021  0.000000e+00
27    AfD - Zweitstimmenanteil: Veränderung 2021-2017  0.000000e+00
7                                      sgbII_pro_kopf  0.000000e+00
15   GRÜNE - Erststimmenanteil: Veränderung 2021-2017  0.000000e+00

xgb_reg.fit(daten_tree, y_tree)

full_pred_boost = xgb_reg.predict(daten_tree)
rmse_full = mean_squared_error(y_tree, full_pred_boost, squared=False)

print(f"RMSE, no validation: {rmse_full:.3f}")
print(f"RMSE, yes validation: {rmse_validation:.3f}")

RMSE, no validation: 0.000
RMSE, yes validation: 6.104

XGBoost ist noch schlechter als der random forest. Außerdem zeigt der Unterschied zwischen den Fehlern mit und ohne Validieren auf frischen Daten einen wichtigen Unterschied zwischen gradient boosting und random forests:
    Es ist im Prinzip unmöglich, random forests überanzupassen, wenn man die einzelnen trees nicht zu tief macht; mehr ist besser. Dasselbe gilt nicht für boosting! Ohne Validieren geht der Fehler des Modells auf Null, weil es sich genau an die Trainingsdaten anpasst.

Was könnte man sonst noch machen? ¶

Hoch